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ABSTRACT 

The expansion of a superbubble is investigated both analytically and numeri- 
cally. Our model implements the thin layer approximation in a vertical profile 
of density as given by an isothermal self-gravitating disk. A precise comparison 
with the results of numerical hydro-dynamics is given. Analogies are drawn 
with the Kompaneets equation that includes the quadratic hyperbolic-secant 
law in the list of the plane-parallel stratified media. An astrophysical applica- 
tion is made to the superbubble connected with the two worms 46.4+5.5 and 
39.7+5.7. The effects of the rotation of the galaxy on the simulated radius and 
on the velocity are introduced. The worms with their strong limb-brightening 
visible on astronomical maps are explained in the framework of image theory. 
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1 INTRODUCTION 



The superbubble (SB) plays a relevant role in astrophysics because (i) it transports ma- 
terial from the galactic plane up to a great galactic height , (ii) it can be a site for th e 
acceleration of cosmic rays, see Higdon fc Lingenfelter ( 2OO5I . 2006 ); Butt fc Bykov ( 2008 ); 



Ferrand &: Marcowithl (120101 ). SBs have been observed in variou s bands: fo r the H I maps, 



sec 



Oey et al 



the X-ray maps, see 



( 12002), for the optical and HI observations, see 



Oev et al 



(2002 



Chul (120081 ): iRodriguez-Qonzalez et all (1201 ll ): I.Taskot et all (1201 lh : of- 



, and for 



ten the image shows a strong limb-brightening, which indicates that the emitting layer is 
thin. The worm is another observed feature that may, or may not, be associated with a 
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wall of an SB. Galactic worms were first identified as irregular, vertical columns of atomic 
gas stretching from the galactic plane; now, simi 



uum a nd infrared maps, see for example 



Koo et al. 



ar structures are found in radio contin- 



f )1992l ): English et all (I2000h : 



BaeketaL 



(120081 ). The models that explain the SB as being due to the combined explosions of su- 
pernova in a cluster of massive stars will now be briefly reviewed. The hydro dynamical 
approximation, with the in clusion of inters t ellar d ensity gradien t s, can produce a blowout 



into the galactic halo, see 



Mac Low et al. 



(1198a ): 



Melioli et al. 



(120091 ) . Expansion in the 



pre sence of m a gnetic fields has been implemente d in various magneto-hydrodynamic codes, 



sec 



Tomisakal (119921 ): iRafikov k Kulsrudl (120001 ) . In semi-analytical calculations, the thin 



layer approxim a tion can be the key to obt a ining the expansion o f the S B: see, for example, 



McCravl (119871 ) 



McCrav k Kafatosl (119871 ): Mac Low k McCravl (I1988h . The thin layer ap 



proximation allows of findin g the equation of motion for an expansion in the framework of 



momentum conservation, see 



Dyson. J. E. and Williams. D. A.I ( 119971 ). when the density of 



the surrounding medium is constant. The case of an expansion in a medium with variable 
density i s more complex and an exponential and a power law vertical profile have been anal- 
ysed, see lZaninettil ( 120101 ) . The exponential and vertical profiles in density do not correspond 
to some physical process of equilibrium. The case of an isothermal self-gravitating disk (ISD) 
is an equilibrium vertical profile which can be coupled with momentum conservation. The 
models cited leave some questions unanswered or only partially answered: 

• Is it possible to calibrate the vertical profile of an ISD? 

• Is it possible to deduce an analytical formula for the temporal evolution of an SB in 
the presence of a vertical profile density as given by an ISD? 

• Is it possible to deduce numerical results for an SB when the expansion starts at a given 
galactic height? 

• What is the influence of galactic rotation on the temporal evolution of an SB? 

• Can we explain the worms as a particular effect using image theory applied to SBs? 

In order to answer these questions, Section [2] reviews the standard equation for momentum 
conservation in an advancing shell, Section [3] introduces a vertical profile in the number of 
particles as given by an ISD which models an aspherical expansion, Section H] applies the 
new law of motion to the SB associated with GW 46.4+5.5, and Section |5] contains detailed 
information on how to build an image of a SB as seen from the equatorial plane as well from 
the poles. 
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2 THE SYMMETRICAL THIN LAYER APPROXIMATION 

The thin layer approximation assumes that all the swept-up gas accumulates infinitely in a 
thin shell just after the shock front. The conservation of radial momentum requires that 



-7rR 3 pR = M 
3 



(i) 



where R and R are the radius and the velocity of the advancing shock, p the density of the 
ambient medium, M n the momentum evaluated at t = tn, R n the initial radius, a nd Rq the 



initial velocity, see 
of motion is 

R = Ro ( 1 + 4 ^( t " t o) 
and the velocity 

R = Ro fl + 4^(t-t ) 



Dvson. J. E. and Williams. D. A.I (119971 ): IPadmanabhanl (I200l[ ). The law 



(2) 



(3) 



3 ASYMMETRICAL LAW OF MOTION 

Given the Cartesian coordinate system (x, y, z), the plane z = will be called the equatorial 
plane, z = Rsin(6 l ), where 9 is the latitude angle which has range [—90° <r» +90°], and R 
is the distance from the origin. The latitude angle is often used in astrophysics to model 
asymmetries in the polar lobes, see the example of the nebula around 77-Carinae (Homuncu- 
lus) shown in Table 1 in 



Smith 



(120021 ) . In our framework, the polar angle of the spherical 
coordinate system is 90 — 9. The vertical number density distribution of galactic HI is usu- 
ally modeled by the following three component behavior as a function of the galactic height 



z which is the distance from the galactic plane in pc: 
n(z) 



nie - Z W + n2e - Z W + n3e -NI/H3 



(4) 



We set the densities in Eq. (J3J) to n x = 0.395 particles cm -3 , n 2 = 0.107 particles cm" 3 , 
n 3 = 0.064 par t icles cm" 3 , and the scale he i ghts t o Hi = 127 pc, H 2 = 318 pc, and H 3 = 403 
pc. ( iLockmanl (119841 ); iDickey fc Lockmanl ( 119901 ); iBisnovatyi-Kogan fc Silichl (Il995l )). This 
distribution of galactic HI is valid in the range 0.4 ^ R ^ R , where Ro = 8.5 kpc and R is 
the distance from the galaxy centre. Here, conversely, we adopt the density profile of a thin 
self-gravitating disk of gas which is characterized by a Maxwellian distribution in velocity 
and distribution which varies only in the z-direction (ISD). The number density distribution 



is 
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Figure 1. Profiles of density versus scale height z: the self-gravitating disk as given by Eq. J5} when h = 90 pc (dashed) and 
the three-component exponential distribution as given by Eq. (J4j (full line). 



n(z) = n sech 2 ( — ) 



(5) 



where iip is the density at z = 0, h is a scaling parameter, and sech is the hyperbolic secant 



( ISpitzerl (119421 ); iRohlfsl ( 119771 ): iBertinl ( 1200(1 ): iPadmanabhanl fl2002h l 

Fig. ([T]) compares the empirical function sum of three exponential disks and the theoret- 
ical function as given by Eq. (jSJ). Assuming that the expansion starts at z = 0, we can write 
z = Rsin(6 l ), and therefore 



n(R, 6) = n sech 



2 Rsin(0). 



2h 



(6) 



where R is the radius of the advancing shell. 

The 3D expansion that starts at the origin of the coordinates will be characterized by 
the following properties. 

• The dependence of the momentary radius of the shell on the latitude angle 9 over the 
range [-90° <R> +90°]. 

• The independence of the momentary radius of the shell from 0, the azimuthal angle in 
the x-y plane, which has a range [0° O 360°]. 



The mass swept, M, along the solid angle A Q between and R is 
M(R,0 
where 



Ail 4 n 

— — m H n I m (R) + -7rR n m H 



I m (R) 



rW(^ffl)dr 



(?) 



(8) 



where Rq is the initial radius and hih the mass of hydrogen. The integral is 



I m (R) = — 4hr 2 (siii(0))- 1 (l +e^ 



-i 



hr 2 
sin(0) 



-8h 2 rln(l +e^)(sin(0))- 2 -8h 3 Li 2 (-e ££ ^)(sin(0))- 3 



(9) 
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where 

Li 2 (z) = E 



oo z n 



n=l 



(10) 



is the dilogarithm, see iHilll ( 118281 ) : iLewinl (Il98l[ ): 



Olver et al. 



The conservation of momentum along the solid angle A Q gives 
M(R,0)R(0) = M(R )R , (11) 

where R(0) is the velocity at R and R is the initial velocity at R = R . Using the previous 
equation, an analytical expression for R(#) along the solid angle can be found, but it is 
complicated, and therefore we omit it. In this differential equation of the first order in R, 
the variables can be separated and integrating term by term gives 



R 



M(r,#)dr = M(R )R (t - t ) , 

Ro 

where t is the time and t the time at R . We therefore have an equation of the type 
1 



(12) 



J-"(R, Ro, h) NL = -R Ro (t — t ) , (13) 

where J-"(R, Ro, 1i)nl has an analytical but complicated form. The case of expansion that 
starts from a given galactic height z, denoted by zqb; which represent the OB associations, 
cannot be solved by Eq. (TT3]) . which is derived for a symmetrical expansion that starts at z = 
0. It is not possible to find R analytically and a numerical method should be implemented. 
In our case, in order to find the root of the nonlinear Eq. (f!3l) . the FORTRAN subroutine 



ZRIDDR from 



Press et al. 



(119921 ) has been used. 



V 



The following two recursive equations are found when momentum conservation is applied: 

Rn+l = Rn + V n At 

— v -(rab) ■ (14) 

where R n , V n , M n are the temporary radius, the velocity, and the total mass, respectively, At 
is the time step, and n is the index. The advancing expansion is computed in a 3D Cartesian 
coordinate system (x, y, z) with the centre of the explosion at (0,0,0). The explosion is better 
visualized in a 3D Cartesian coordinate system (X, Y, Z) in which the galactic plane is given 
by Z = 0. The following translation, Tob, relates the two Cartesian coordinate systems. 
X = x 



OB 



Y = y 

Z = Z + Zqb 



(15) 
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where zqb is the distance in parsecs of the OB associations from the galactic plane. 

The physical units have not yet been specified: parsecs for length and 10 7 yr for time are 
perhaps an acceptable astrophysical choice. With these units, the initial velocity V = Ro 
is expressed in units of pc/(10 7 yr) and should be converted into km/s; this means that 
Vq = 10.207Vi where Vi is the initial velocity expressed in km/s. 



Analytical results can also be obtained solving the Kompaneets equation, see 



Kompaneets 



( 1l960l ). for the motion of a shock wave in different plane-parallel stratified media such as 
exponential, power-law type, and a quadratic hyperbolic-secant, see synoptic Table 4 in 



Olanol ((20091). 



4 ASTROPHYSICAL APPLICATIONS 



A use ful formula that allows setting up the initial conditions is formula (10.38) in iMcCray 
( 119871 ) which models the SB evolution in the energy conserving phase when the number 
density is constant, 

R= 111.56 pc(5^)*, (16) 
no 

where is the time expressed in units of 10 7 yr, E51 is the energy expressed in units of 
10 51 erg, n is the number density expressed in particles cm -3 (density p = n m, where 
m = 1.4mn) and N * is the number of SN explosions in 5.0 • 10 7 yr and therefore is a rate, 



sec 



McCrayl (119871 ). This formula, deduced in spherical coordinates, can be used only when 



the density is constant, as an example for the first 20 pc, the variation of the density with 
galactic height is ~ 2%. In the following, we will assume that the bursting phase ends at 
t = t^o (the bursting time is expressed in units of 10 7 yr) when Nsn SN are exploded 

N SN = N*^f . (17) 
5 • 10 K J 

The velocity of the SB in the energy conserving phase when the number density is 
constant is 



6.567 s M51Ml k 

v(to= t y<> v • (18) 

Eq. (TIB]) and (IT51) can be used to deduce R and and V once N* and t^p are given. We 
continue giving some information about the astrophysical target of our simulation, analysing 
the analytical solution at zqb =0 and the numerical solution zqb 7^ 0, making a comparison 
with the hydro code and analysing the effects of the galactic rotation on the obtained results. 
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Table 1. Data of the SB associated with GW 46.4+5.5. 



Dimensions (pc ) 
Expansion velocity (km s 
Age (10 7 yr) 
zob (pc) 
Total energy (10 51 erg) 



- 1 ) 



345 x 540 
15 
0.5 
100 
15 




-200 0„ , * 200 400 

Y (pc) 

Figure 2. Section of the SB GW 46.4+5 in the Y-Z;X=0 plane when the explosion starts at zob = Opc. The analytical 
results arc shown as full line and the numerical results as crosses. The code parameters of the solution of Eq. (1 1 3 P as well of 
the numerical couple Q3} are h = 90 pc, t 7 = 0.45 , t 7 , = 0.0045, r = 49.12, V = 641.7kms- 1 , N SN = 93 and N*=103000. 



4.1 The simulated object 



When the two worms 46.4+5.5 and 39.7+5.7 are carefully studied (jKim Sz Kool (120001 )) it is 
possible to conclude that they belong to a single SB. The parameters of this single SB are 
given in Table [U 



4.2 Analytical versus numerical solutions 

The analytical model as given by the solution of the nonlinear Eq. ( 113]) can be used only in 
the case zqb=0. As an example, we give a model of the SB associated with GW 46.4+5.5, 
see Fig. [2j The numerical solution as given by the recursive relationship (fT4"]) can be found 
adopting the same input data of the analytical solution, see crosses in Fig. (F5]); the numerical 
solution agrees with the analytical solution within 0.65%. We are now ready to present the 
numerical evolution of the SB associated with GW 46.4+5 when zob = 100pc, see Fig. [3] 
and Table [2j 

Table 2. Simulated data of the SB associated with GW 46.4+5.5. 

Size (pc 2 ) 454 ■ 521 

Averaged expansion velocity (km s -1 ) 19.3 
Age (10 7 yr) 0.45 
zob (pc) 100 
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Y (pc) 



Figure 3. Section of the SB GW 46.4+5 in the Y-Z;X=0 plane when the explosion starts at zqb = lOOpc. The code 
parameters for the numerical couple (1141 are h = 90 pc, t 7 = 0.45 , t 7>Q = 0.00045, r = 24.43, V = 3191kms-\ N SN = 180 
and N* =2000000. The explosion site is represented by a cross. 



4.3 Numerical solution and hydro code 

The level of confidence in our results ca n be given by a comp arison with numerical hydro 
dynamics calculations (see, for example 
bution they adopted, see equation (1) in 



Mac Low et a 



Mac Low et al 



(1989)). The vertical de nsity distri- 



(Il989l ) and equation (5) in 



Tomisaka fc Ikeuchi 



(119861 ) . has the following dependence on z, the distance from the galactic plane in parsec: 



nhydro = n d <^ 6exp[- 



Vjz) 



a 



+ (1 - 9)exp 



ic 



Vp(z) 



with the gravitational potential 



V p (z) = 68.61n 



1 + 0.9565 sinh 2 0.758 



z 



(kms 



-1\2 



(19) 



(20) 



Here, n^ = 1 particles cm -3 , = 0.22, a\c = 14.4km s _1 , o"c = 7.1km s _1 , and z = 124 pc. 
Fig. 01]) compares the hydro number density as given by Eq. [19] and the theoretical function 
as given by Eq. ([3]). The difference in the density profiles from hydrosimulations and from 
the simple model adopted in this paper are due to the fact that t he density at z = is 



assumed to be n h . 



ydro = 1 particles cm 3 in the hydrosimulations, see 
In our model conversely at z = we have n = 0. 5 66 pa rticles cm -3 as in 



Mac Low et al. 



A typical run of ZEUS (see 



Mac Low et al. 



Lockmanl (119841 ). 



1989). 



(]1989D ). a two-dimensional hydrodynamic 



code, is done for tj = 0.45 and supernova luminosity of 1.6 ■ 10 38 erg s _1 when zob=100 pc. 
In order to make a comparison with our code, we adopt the same time and we search the 
parameters which produce similar results, see Fig. [5] The percentage of reliability of our 
code can also be introduced, 
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Figure 4. Profiles of density versus scale height z: the self-gravitating disk as in Eq. {5} when h = 90 pc (dashed line) and 
the hydro number density as given by Eq. 1191 (full line) . 




Y (pc) 



Figure 5. Section of the SB in the X-Z plane when the explosion starts at zob = lOOpc (empty stars). The code parameters 
are h = 90 pc, t 7 = 0.45 , t 7 ,o = 0.00045, r = 24.43, V = 3191 kms" 1 , N SN = 180 and N*=2000000. The points represented 
by the small crosses indicate the inner section from Fig. 3a of MacLow et al. 1989. The explosion site is represented by a big 
cross. 



e=(1 _ KRbydro Rnu m )| ) lQ0) (21) 
-^hydro 

where Rhydro is the radius given by the hydro-dynamics, and R num s the radius obtained from 
our simulation. Table [3] shows our numerical radii in the upward, downward, and equatorial 
directions, and the efficiency as given by formula ( 12T1) . 



Table 3. Code reliability. 





Rup(pc) 


Rdown(P c ) 


Req(pc) 


Rhydro (ZEUS) 


330 


176 


198 


Rnum (our code) 


302 


203 


217 


efficiency(%) 


91 


84 


90 
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Figure 6. Section of the SB GW 46.4+5 in the Y-X,Z=0 plane when the explosion starts at zob = lOOpc. Full points represent 
the circular section, the empty stars the rotation-distorted section. Parameters as in Fig. [3] and Ro[pc]=8500. 



4.4 Galactic rotation 



The influence of the Galactic rotation on the results can 
law of the Galactic rotation as given by 

-Ro[pck o.382i 



Wouterloot et al 



) e be o btained by introducing the 
(1199oL 



Vr(R c 



220 



8500 



2 km sec 1 , 



(22) 



where Ro is the radial distance from the centre of the Galaxy in parsecs. The original circular 
shape of the superbubble at a given value of z transforms to an ellipse by the following 
transformation, T r , 



T r <^ 



x/ = x + 0.264y t 

Y' = Y 
z/ = z, 



(23) 



Zaninetti! (120041 ) . see Fig. El In the 



where y is in parsecs and t is expressed in units of 10 7 yr 
same way, the effect of the shear velocity as a function of the distance y from the centre of 
the expansion, V sm ft(y), can be easily obtained by performing a Taylor expansion of Eq. 



y 



V s hift(y) = 84.04^- Km/s 
K 



(24) 



This is the the shear velocity as a function of y. Fig. [7] shows the rotation-distorted map of 
position and velocity. 

The effect of rotation seems to be small and in order to see a more relevant correction 
we should increase the elapsed time, see Fig. |HJ 
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Figure 7. The stars represent the rotation-distorted section when Ro[pc]=8500. The arrows show the velocity field of the 
expansion as modified by the shear velocity. The direction of rotation of the Galaxy (big arrow) and the scale of the velocity 
are shown. 




Figure 8. Section of the SB GW 46.4+5 in the Y-X,Z=0 plane when the explosion starts at zob = lOOpc. Full points represent 
the circular section and the empty stars the rotation-distorted section. Parameters as in Fig.|3]except t-? = 2 and Ro[pc]=8500. 



5 THE IMAGE 



he tr an sfer equation in th e 



fll99lh or 



Hjellming. R. M 



prese nce of emission only, see for example iRybicki fc Lightman 



(119881 ). is 



^ = -KO-v + ji/C 
ds 



(25) 



where \ v is the specific intensity, s is the line of sight, ] v the emission coefficient, k u a mass 
absorption coefficient, ( the mass density at position s, and the index v denotes the relevant 
frequency of emission. The solution to Eq. (I25p is 



Ut v ) 



where t v is the optical depth at frequency v 
dr u = k u (ds . 



(26) 



(27) 



12 L. Zaninetti 

We now continue analysing the case of an optically thin layer in which r v is very small ( or 
k u very small ) and the density ( is replaced with our number density C(s) of particles. One 
case is taken into account: the emissivity is proportional to the number density 

U = KG(s) , (28) 

where K is a constant. This can be the case for synchrotron emission and a simple model 
for the acceleration of the electrons can be found in Appendix |A] We select as an exam- 
ple the [S II] con t inuum of the synchrotron superbubble in the irregular galaxy IC10, see 



Lozinskaya et al. 



( 120081 ). and the X-ray emission below 2 keV around the OB assoc iation 



LH9 in the H II complex Nil in the Large Magellanic Cloud, see 
intensity at a given frequency is 

oc lz/ , 



Maddox et al. 



fl2009h . The 



(29) 



where 1 is the length of the radiating region along the line of sight. The source of synchrotron 

luminosity is assumed here to be the rate of kinetic energy, L m , 
1 



L, 



-pAV d 



(30) 



where A is the considered area, see formula (A28) in lDe Young! (120021 ). In the case of formula 
(TT6]) A = 47rR 2 , which means 

L m = ip47rR 2 V 3 , (31) 

where R is the instantaneous radius of the SB and p is the density in the advancing layer 
in which the synchrotron emission takes place. The astrophysical version of the the rate of 
kinetic energy, 



-,29 D 2 3 er § s 



1.39 x 10 zy rit RSVj 



(32) 



where ni is the number density expressed in units of 1 p , Ri is the radius in parsecs, 
and Vi is the velocity in km/s. The spectral luminosity, L u , at a given frequency v is 

L u = 4ttD 2 S i/ , (33) 

with 



v 



(34) 



Sjy — Sn( — 

where So is the flux observed at the frequency vq and D is the distance. The total observed 
synchrotron luminosity, L tot , is 



L 



tot 



(35) 
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Table 4. Table of the values of f c when u min = 10 7 Hz, v max = 10 18 Hz and = -0.7. 



band 


A (A) 


FWHM (A) 




fc 


U 


3650 


700 


6.86 


xlO -3 


B 


4400 


1000 


7.70 


xlO" 3 


V 


5500 


900 


5.17 


xl0~ 3 


Ha 


6563 


100 


0.56 


xl0~ 3 


[SII] continuum 


7040 


210 


0.92 


xl0~ 3 




Figure 9. The two circles (sections of spheres) which include the region with constant density are represented through a full 
line. The observer is situated along the x direction, three lines of sight are indicated, and the relativistic electrons have radius 
r in the region a < r < b. 

where z/ m ; n and v m&K are the minimum and maximum frequencies observed. The total ob- 
served luminosity can be expressed as 



L 



tot 



(36) 



where e is a constant of conversion from the mechanical luminosity to the total observed 
luminosity in synchrotron emission. 

The fraction of the total luminosity deposited in a band f c is 
v ■ /3+1 - v 



f, 



P+l _ (3+1 
'c,mm f c.max 



^min u max 

where z/ Cjm ; n and z/ C;max are the minimum and maximum frequency of the band. Table [5] shows 
some values of f c for the most important optical bands. An analytical solution for the radial 
cut of intensity of emission can be found in the equatorial plane zqb = 0. We assume that 
the number density of relativistic electrons C is constant and in particular rises from at 
r = a to a maxi mum value C m , remains constant up to r = b, and then falls again to 0, see 
Zaninettil (120091 ). This geometrical description is shown in Fig. [9j The length of sight, when 
the observer is situated at the infinity of the x-axis, is the locus parallel to the x-axis which 
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pc from center 

Figure 10. Cut of the mathematical intensity / of the ring model, Eq. I|38|l , in the equatorial plane (full line) of the SB 
associated with GW 46.4+5.5 . The x- and y-axes are in pc, a = 172.5 pc, b = 186.9 pc and j^Zq| = 5. 

crosses the position y in a Cartesian x — y plane and terminates at the external circle of 
radius b. When the number density of the relativistic electrons C m is constant between the 
two spheres of radii a and b, the intensity of radiation is 

I 0a = C m x 2 x (^/b 2 - y 2 - ^a 2 - y 2 ) ;0^y<a 

; a ^ y < b . (38) 



iab = c m x 2 x (v b2 - y 



The ratio between the theoretical intensity at the maximum ((y = a)) and at the minimum 
(y = 0) is given by 
Vb 2 - a 2 



i(y = a) 



(39) 



I( y = 0) b - a 

A cut in the theoretical intensity of the SB associated with GW 46.4+5.5 is shown 



in Fig. [TOj Similar analytical results for the intensity of the cuts in t he Ha of p 



nebulae and in the radio of supernova remnants have been found by iGrav et al. 



anetary 

(bom 



compare their Figure 5 with our Figure and by lOpsenica fc Arbutinal (120111 ) , see their 
Figure 1. A simulated image of the complex shape of a SB is composed by combining the 
intensities which characterize different points of the advancing shell. For an optically thin 
medium, the transfer equation provides the emissivity to be multiplied with the distance on 
the line of sight, 1. This length in SB depends on the orientation of the observer but for the 
sake of clarity the observer is at infinity and sees the SB from the equatorial plane zq-b =0 
or from one of the two poles. We now outline the numerical algorithm which allows us to 
build the complex image of an SB. 

• An empty (value=0) memory grid .M(i,j,k) which contains NDIM 3 pixels is considered 

• We first generate an internal 3D surface by rotating the section of 180° around the 
polar direction and a second external surface at a fixed distance AR from the first surface. 
As an example, we fixed AR = 0.03R max , where R max is the maximum radius of expansion. 
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pc from center 



Figure 11. Map of the theoretical intensity of SB GW 46.4+5. Physical parameters as in Fig. [2] The three Euler angles 
characterizing the orientation are $=90 °, ©=90 ° and \I/=90 °. This combination of Euler angles corresponds to the rotated 
image with the polar axis along the z-axis. In this map It r = Imax/2. 



The points on the memory grid which lie between the internal and the external surfaces 
are memorized on A^(i,j,k) with a variable integer number according to formula (T3TTT) and 
density p proportional to the swept mass. 

• Each point of ,M(i,j,k) has spatial coordinates x, y, z which can be represented by the 
following 1x3 matrix, A, 



x 



A 



The orientation of the object is cha racterized by the 



a total 3x3 rotation matrix, E, see 
by the following 1x3 matrix, B, 

B = E • A . 



Goldstein et al. 



(40) 



uler angles ($, 0, \I/) and therefore by 



( 120021 ). The matrix point is represented 



(41) 



• The intensity map is obtained by summing the points of the rotated images along a 
particular direction. 

• The effect of the insertion of a threshold intensity, I tr , given by the observational tech- 
niques, is now analysed. The threshold intensity ca n be parametrize d by I maX ; the maximum 



value of intensity which characterizes the map, see 



Zaninetti 



feojj ) 



An ideal image of the intensity of SB GW 46.4+5 is shown in Fig. [H} and Fig. [12] shows 
two cuts through the centre of the SB. 

We can also build the theoretical image as seen from one of the two poles, see Fig. | 
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R (fpc) 



Figure 12. Two cuts of the model intensity across the centre of the SB: equatorial cut (full line) and polar cut (dotted line). 
Parameters as in Fig. Illl and Itr = 0. 




-400 -200 200 400 

pc from center 

Figure 13. Map of the theoretical intensity of SB GW 46.4+5 as seen from the polar direction at infinity. Physical parameters 
as in Fig. [2] The three Euler angles characterizing the orientation are $=0°, ©=0° and ^=0°. In this map Itr — Imax/2. 



6 CONCLUSIONS 

Law of motion The temporal evolution of an SB in a medium with constant density is 
characterized by a spherical symmetry. The presence of a thin self-gravitating disk of gas 
which is characterized by a Maxwellian distribution in velocity and distribution of density 
which varies only in the z-direction produces an axial symmetry in the temporal evolution 
of an SB. The resulting Eq. ffTB"]) has an analytical form which can be solved numerically 
when z b =0. The case of z G b 7^ can be attached solving two recursive equations, see (JHj). 
These complex shapes can also be modeled by the Kompaneets equation when a quadratic 
hyperbolic-secant is adopted. Both the approximation here used and the Kompaneets equa- 
tio n can model the shape of the SB without invoking the collision of two wind-blown SBs, 



sec 



Ntormousi et al. 



(120111 ). The code developed here runs in less than a minute on a LINUX 
-2.66GHz processor and can be an alternative to a purely numerical 2D or 3D model. 
Synchrotron Emission 

Here we assumed that the conversion from the flux of kinetic energy into non thermal 
luminosity is regulated by a constant. The exact value of this constant can be deduced 
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once the number density, radius, thickness, and velocity of the advancing shell are provided 
by observational astronomy. The values of the fraction of flux deposited in the various 
astronomical bands, f c , are given in Table [5j 
Images 

The emissivity in the thin advancing layer is assumed to be proportional to the rate of 
kinetic energy, see Eq. (130]) . where the density is assumed to be proportional to the swept 
material. This assumption explains the strong limb-brightening visible on the astronomical 
maps of SBs and allows associating the observed worms to the non detectable SBs. As an 
example, if the threshold of the observable flux at a given wavelength is bigger than the 
emitted flux at the centre of an SB, we will detect only the external regions, the so called 
worms, see Fig. [T2l 
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APPENDIX A: HOW TO ACCELERATE ELECTRONS 

An electron which loses its energy due to synchrotron radiation has a lifetime of 
E 

Tr ~ ~ 500E -1 H- 2 sec , (Al) 
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where 



power ( jLang 



is the en ergy in ergs, H the magnetic field in Gauss, and P r is the to tal radiated 



19991 . Eq. 1.157). The energy is connected to the critical frequency (jLang 



1999, 



Eq. 1.154), by 

u c = 6.266 x 10 18 HE 2 Hz . 

The lifetime for synchrotron losses is 
1 



syn 



39660 



HVH^ 



yr 



(A2) 



(A3) 



Following ( jFermi 



1949 



19541 ) . the gain in energy in a continuous form for a particle which 

(A4) 



spirals around a line of force is proportional to its energy, E, 
dE _ E 
dt r n 

where rn is the typical time-scale, 
1 4 u 2 c 

- = |( ] r)(r-) > 

Tn 3 c 2 L n 



(A5) 



where u is the velocity of the accelerating cloud belonging the adva ncing shell of the SB, c 



is the velocity of light, and L n is the mean free path between clouds (jLang 



19991 . Eq. 4.439). 



The mean free path between the accelerating clouds in the Fermi II mechanism can be found 
from the following inequality in time: 

Til < Tsync , (A6) 

which corresponds to the following inequality for the mean free path between scatterers 

(A7) 



1.72 10 5 u 2 

L < ; PC 



HVHz/c 2 

The mean free path length for an SB at [S II] continuum, which means 7040 A or 704 nm 
or 4.258 10 14 Hz, gives 



L < 



2.94 10" 6 ut' 



pc 



(A8) 



where ui is the velocity of the accelerating cloud expressed in km/s, and H5 = is the mag- 
netic field expressed in units of 10 -5 Gauss. When this inequality is verified, the direct 
conversion of the rate of kinetic energy into radiation can be adopted. Recall that the Fermi 
II mechanism produces an inverse power law spectrum in the energy of the type N(E) oc E~ 7 



or an inverse p ower law in the observed frequencies N(z/) oc with /3 = —- L n- ( Lang 



1999 



Zaninei 


tti 


2011) 


by 


Bell 


( 


L978a) 



2011 ). Th e str o ng sho ck accelerating mechanism named Fermi I was introduced 



Belli (jl978al ) and iBelll ( jl978bl ). The energy gain relative to a particle that is crossing the 



shock is: 
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(A9) 



2 u. , c . 

3V /It 
c Li 



(A10) 



dE _ E 

dF ~ Ti ' 

where T\ is the typical time-scale, 
1 

n 

where u is the velocity of the shock of the SB, c is the velocity of light, and Lj is the mean 
free path between scatterers. This process produces an energy spectrum of electrons of the 
type: 

N(E)dE oc E" 2 dE 



(All) 



sec 



Longairl (I1994j ). The two mechanisms , Fermi I and Fermi II, produce the same results 
when 



Ln 
Li 



u 



(A12) 



